/*
**  File: AquaFuncts_biosym.c
**
**  This file is part of the AQUA program suite, which is being
**  developed in the NMR Spectroscopy Department, Bijvoet Center for 
**  Biomolecular Research, Utrecht. 
**
**
**  AUTHOR(S) : Ton Rullmann
**  VERSION   : 1.1
**  DATE      : Feb 17, 1997
**
**  This file contains routines for handling Biosym files.
**
**  Most header definitions are present in AquaFuncts_io.c
*/

/****************  function prototypes  ****************/

#include "AquaFuncts_biosym.h"

static Boolean ReadHeader_biosym( FILE *fp, char *header, const size_t size, 
                                  const Boolean first_frame );
static Status ReadAtomAndCoordData_biosym( FILE *fp, char *atnam, char *resnam, 
        int *atnum, int *resnum, float *px, float *py, float *pz, 
	const size_t size, const int count );
static Status ReadAtomData_biosym( FILE *fp, char *atnam, char *resnam, 
                                   int *atnum, int *resnum, const int size );
static Status ReadCoordData_biosym( FILE *fp, float *x, float *y, float *z );
static void MDFStorePseudo( char *line, Psatom *psatom );
static void MDFStoreAtom( char *line, Psatom *psatom );

/* BEGIN_EXTH */

/****************  ReadFil_cv_biosym  ****************/
/*
**  Read a Biosym car file in the ReadFil_cv way.
**
*/

Status ReadFil_cv_biosym( FILE *fp, const Filtyp type, const int atom_count,
                          int *fil_struc_num, int *nerr )
{
    int      count = 0;
    char     *a, *r;
    int      *ia, *ir;
    float    *px, *py, *pz;
    Status   stat;

/* allocate memory */

    Alloc_cv( atom_count );
        
/* read the file header - for all structures */
/* but test the first line only for the first structure */

/*    if ( !*fil_struc_num )   */

        ReadHeader_biosym( fp, header, size_HEADER, (*fil_struc_num? false:true) );
    
/* read the predetermined number of atoms */

    (*fil_struc_num)++;
    a = atonam;
    r = resnam;
    ia = atono;
    ir = resno;
    px = x;
    py = y;
    pz = z;

    while ( ++count <= atom_count )
    {
        stat = ReadAtomAndCoordData_biosym( fp, a, r, ia++, ir++, 
                                            px++, py++, pz++, size_ID, count );
        if ( stat != status_ok )
	{
	    if ( stat == eof_err && count == 1 )
		return( end_of_data );
	    else
                goto err;
	}
        a += size_ID;
        r += size_ID;

        if ( FlushLine( fp ) < 0 )
	{
	    stat = eof_err;
            goto err;
	}
    }
    return( status_ok );

err:
    *nerr = count-1;
    puts( " * Error or unexpected end of file" );
    return( stat );
}

/****************  AtomCount_biosym  ****************/
/*
**  Count the number of atoms per structure in a Biosym car or arc file.
*/

int AtomCount_biosym( FILE *fp )
{
    Card     line;
    int      count = 0, endcnt = 0;
    
    if ( !ReadHeader_biosym( fp, header, size_HEADER, true ) )
	return( 0 );
    while ( ReadLine( line, size_LINE, fp ) )
    {
	if ( EQn( line, "PBC", 3 ) || EQn( line, "HELIX", 5 ) )
	    continue;
	if ( EQn( line, "end", 3 ) )
	{
	    if ( ++endcnt == 2 )
		break;           /* two consecutive 'end' lines signal end of frame */
	    else
		continue;
	}
	endcnt = 0;              /* reset 'end' count */
	count++;
    }
    return( count );
}

/****************  NoeCount_biosym  ****************/
/*
**  Determine the number of Noe restraints from a Biosym data file
*/

int NoeCount_biosym( FILE *fp )
{
    char   *fct_name = "NoeCount_biosym";

#ifdef FPOSON
    int    count = 0;
    Card   line;
    int    c;
    fpos_t *file_pos = NULL;
    
    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( c != '#' )
        {
            FlushLine( fp );
            continue;
        }
        if ( ReadLine( line, size_LINE, fp ) )
            if ( EQn( line, "NOE_distance", 12 ) )
            {
                fgetpos( fp, file_pos );
                break;
            }
    }
    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( c == '#' )
            break;
        else
            if ( ReadLine( line, size_LINE, fp ) )
                if ( LenChars( line, size_LINE ) )
                    count++;
    }
    if ( !file_pos )
	FATAL( "could not determine position within file" );

    fsetpos( fp, file_pos );  /* doesn't seem to work........ */

    return( count );
#else
    FATAL( "Sorry, this routine could not be compiled" );
    return(1); /* This statement will never be executed */
#endif
}

/****************  ReadNoes_biosym  ****************/
/*
**  Read Noe restraints from a Biosym data file. (cf. ReadRestraintFile)
**  The restraint atom names are stored in the external field.
**  The residue ID is completed ** by calling MakeUniqResid **.
**
**  Biosym format:
**
**  #NOE_distance
**  1:PRO_2:HA         1:ASP_3:HN          2.000  4.200  4.200 10.00 10.00 1000.000 0.00
*/

int ReadNoes_biosym( FILE *fp, DRestraint *noes, const int count, Card ident )
{
    int   noe_read = 0;
    int   c, atnum;
    Atom  *at;
    Card  line;

/* position the file */

    rewind( fp );
    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( c != '#' )
        {
            FlushLine( fp );
            continue;
        }
        if ( ReadLine( line, size_LINE, fp ) )
            if ( EQn( line, "NOE_distance", 12 ) )
                break;
    }

/* set ident to empty string */

    *ident = '\0';
    
/* process upto first comment line */

    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( c == '#' )
            break;
        else
        {
            if ( ++noe_read > count )
            {
                printf( " * More input than the expected %i restraints - cannot handle",
                        count );
                exit( 1 );
            }
            at = &noes->atom_A;
            ReadAtomData_biosym( fp, at->atom_ext, at->res_nam, &atnum, 
                                 &at->res_num_orig, size_ID );
            MakeUniqResid( at );
            at = &noes->atom_B;
            ReadAtomData_biosym( fp, at->atom_ext, at->res_nam, &atnum, 
                                 &at->res_num_orig, size_ID );
            MakeUniqResid( at );

            fscanf( fp, " %f %f", &noes->lo_bnd, &noes->up_bnd );
	    noes->stat = bnd_uplo;

            FlushLine( fp );
/*
            WriteAtom( stdout, &noes->atom_A, use_ext );
            WriteAtom( stdout, &noes->atom_B, use_ext );
            printf("%f %f\n",noes->lo_bnd,noes->up_bnd);
*/
            noes++;
        }
    }
    return( noe_read );
}

/****************  ReadMDFPseudo_biosym  ****************/
/*
**  Process the pseudo atom definitions from a biosym mdf file.
*/

#define max_PSAT 500

Psatom *ReadMDFPseudo_biosym( FILE *lib_fp, int *pseudo_count )
{
    Psatom  *psatoms, *pslist;
    int     num_psat = 0;
    Card    line;

    psatoms = CREATE_NEW( Psatom, max_PSAT );
    if ( !psatoms )
    {
	puts( " * Out of memory" );
	exit( 1 );
    }
    InitPsatomS( psatoms, max_PSAT );
    pslist = psatoms;
    
    while ( ReadLine( line, size_LINE, lib_fp ) )         /* input line; eof test */
        if ( EQn( "@list", line, 5 ) )                    /* start of pseudo ? */
        {
            if ( ++num_psat > max_PSAT )
            {
                printf( "\n * Out of array bound, change max_PSAT\n" );
                exit( 1 );
            }
            MDFStorePseudo( line, psatoms );              /* store pseudo */
            while ( ReadLine( line, size_LINE, lib_fp ) ) /* input line; eof test */
            {
                if ( !IsEmptyStr( line ) )                /* blank line as separator */
                    MDFStoreAtom( line, psatoms );
                else
                    break;                                /* stop storing */
            }
            printf( "\nnum_psat %i ", num_psat );
            printf( "number of constituent atoms: %i\n", psatoms->num_form );
            psatoms++;
        }
    *pseudo_count = num_psat;
    return( pslist );
}

#undef max_PSAT

/****************  WriteRestraints_biosym  ****************/
/*
**  Write restraints in Biosym restraint format. 
**  Biosym format:
**
**  #NOE_distance
**  1:PRO_2:HA         1:ASP_3:HN          2.000  4.200  4.200 10.00 10.00 1000.000 0.00
*/

#define TAB 5

void WriteRestraints_biosym( FILE *fp, DRestraint *restrs, int DR_count, 
                             int skip_count, Card ident, char *DR_typ )
{
    char    *fct_name = "WriteRestraints_biosym";
    size_t  l, tot;
    float   kL = 1.0F, kU = 1.0F, maxForce = 1000.0F;

/*
    kL = GetOptionRVal( opt_lowfc );
    kU = GetOptionRVal( opt_uppfc );
    maxForce = GetOptionRVal( opt_maxf );
*/
    fprintf( fp, "!BIOSYM restraint 1\n" );
    fprintf( fp, "!%s\n\n", ident );
    if ( EQ( DR_typ, "noe" ) )
	fprintf( fp, "#NOE_distance\n" );
    else
	fprintf( fp, "#distance\n" );
    while ( DR_count-- )
    {
	if ( !restrs->keep ) 
	{
	    restrs++;
	    continue;
	}
	if ( ( l = LenChars( restrs->atom_A.chain, size_ID ) ) )
	{
	    FPrintChars( fp, restrs->atom_A.chain, l );
	    tot = l;
	}
	else
	{
	    fputc( '1', fp );
	    tot = 1;
	}
	fputc( ':', fp );
	l = LenChars( restrs->atom_A.res_nam, size_ID );
	FPrintChars( fp, restrs->atom_A.res_nam, l );
	fputc( '_', fp );
	tot += ( 2 + l + (size_t)fprintf( fp, "%i", restrs->atom_A.res_num_orig ) );
	if ( restrs->atom_A.res_insert != ' ' )
	{
	    fprintf( fp, "%c", restrs->atom_A.res_insert );
	    tot++;
	}
	fputc( ':', fp );
	l = LenChars( restrs->atom_A.atom_ext, size_ID );
	FPrintChars( fp, restrs->atom_A.atom_ext, l );
	tot += ( l + 1 );
	l = TAB - ( tot % TAB );
	while ( l-- )
	    fputc( ' ', fp );

	if ( ( l = LenChars( restrs->atom_B.chain, size_ID ) ) )
	{
	    FPrintChars( fp, restrs->atom_B.chain, l );
	    tot = l;
	}
	else
	{
	    fputc( '1', fp );
	    tot = 1;
	}
	fputc( ':', fp );
	l = LenChars( restrs->atom_B.res_nam, size_ID );
	FPrintChars( fp, restrs->atom_B.res_nam, l );
	fputc( '_', fp );
	tot += ( 2 + l + (size_t)fprintf( fp, "%i", restrs->atom_B.res_num_orig ) );
	if ( restrs->atom_B.res_insert != ' ' )
	{
	    fprintf( fp, "%c", restrs->atom_B.res_insert );
	    tot++;
	}
	fputc( ':', fp );
	l = LenChars( restrs->atom_B.atom_ext, size_ID );
	FPrintChars( fp, restrs->atom_B.atom_ext, l );
	tot += ( l + 1 );
	l = TAB - ( tot % TAB );
	while ( l-- )
	    fputc( ' ', fp );
	
	switch( restrs->stat )
	{
	    case bnd_upper:
		fprintf( fp, "   %8.3f %8.3f", -1.0, restrs->up_bnd );
		if ( EQ( DR_typ, "noe" ) )
		    fprintf( fp, " %8.3f", restrs->up_bnd );
		break;
	    case bnd_lower:
		fprintf( fp, "   %8.3f %8.3f", restrs->lo_bnd, 9999.0F );
		if ( EQ( DR_typ, "noe" ) )
		    fprintf( fp, " %8.3f", 9999.0F );
		break;
	    case bnd_uld:
	    case bnd_uplo:
		fprintf( fp, "   %8.3f %8.3f", restrs->lo_bnd, restrs->up_bnd );
		if ( EQ( DR_typ, "noe" ) )
		    fprintf( fp, " %8.3f", restrs->up_bnd );
		break;
	    default:
		FATAL( "wrong data type" );
	}
	fprintf( fp, " %8.3f %8.3f %10.3f\n", kL, kU, maxForce );

	restrs++;
    }
}

#undef TAB

/* END_EXTH ----------------------- only private functions follow */

/****************  ReadHeader_biosym  ****************/
/*
**  Read header information from a Biosym car or arc file.
**  Format assumed: 
**  - first line should begin with: "!BIOSYM archive" followed by an integer
**  - text header copied from line immediately before line with "!DATE"
**    (only the first MAXHEADER of such lines are stored) 
**  - data start after line beginning with: !DATE
**    (but PBC or HELIX records may be there; these are not read by this
**    routine)
**
**  After returning the file is positioned at the beginning of a new line.
*/

static Boolean ReadHeader_biosym( FILE *fp, char *header, const size_t size, 
                                  const Boolean first_frame )
{
    Card  line, oldline;
    char  *h = header;
    int   numh = 0, dum;

    if ( first_frame )
    {
	if ( fscanf( fp, "!BIOSYM archive %i", &dum ) < 1 )
	{
	    puts( " * This is not a BIOSYM archive file" );
	    return( false );
	}
    }
    
    while ( ReadLine( line, size_LINE, fp ) )
    {
        if ( EQn( line, "!DATE", 5 ) )
        {
            if ( numh++ < MAXHEADER && alloc_cv )
            {
                strncpy( h, oldline, size );
                h += size;
            }
            break;
        }
	strcpy( oldline, line );
    }
    return( true );
}

/****************  ReadAtomAndCoordData_biosym  ****************/
/*
**  Read atom and coordinate information from one Biosym car file line.
**  Format assumed:
**  N      -32.011642456    1.557512760    6.453227043 VALn   1  n3 N   -0.5000    1
**
**  After reading the file is **not** flushed to the next line.
**  Returns status flag.
*/

static Status ReadAtomAndCoordData_biosym( FILE *fp, char *atnam, char *resnam, 
       int *atnum, int *resnum, float *px, float *py, float *pz, 
       const size_t size, const int count )
{
    Card        line;
    Status      stat;
    static int  endcnt = 0;
/*
    if ( !ReadWord( line, size_LINE, true, fp ) )
	return( eof_err );
*/
    while ( ReadWord( line, size_LINE, true, fp ) )
    {						/* input first word */
	if ( EQn( line, "PBC", 3 ) || EQn( line, "HELIX", 5 ) )
	{
	    FlushLine( fp );
	    continue;
	}
	if ( EQn( line, "end", 3 ) )
	{
	    if ( ++endcnt == 2 )
		break;           /* two consecutive 'end' lines signal end of frame */
	    else
	    {
		FlushLine( fp );
		continue;        /* THERE SHOULD BE A PROVISION TO ASSIGN A NEW CHAIN */
	    }
	}
	endcnt = 0;              /* reset 'end' count */



	CopyChars( atnam, line, size );         /* read atom name */
    
	stat = ReadCoordData_biosym( fp, px, py, pz );
						/* read coordinates */
	if ( stat != status_ok )
	    goto err;
    
	if ( !ReadWord( line, size_LINE, true, fp ) )
	{
	    stat = eof_err;
	    goto err;
	}					/* input next word */
    
	CopyChars( resnam, line, size );        /* read residue name */
	
	if ( fscanf( fp, "%i", resnum ) < 1 )
	{
	    puts( "Error in reading residue number from Biosym file" );
	    stat = read_err;
	    goto err;
	}                                       /* read residue number */

	*atnum = count;                         /* set atom number */

	return( status_ok );
    }
    return( eof_err );

err:
    printf( " * Last read: %s\n", line );
    return( stat );
}

/****************  ReadAtomData_biosym  ****************/
/*
**  Read atom information from one Biosym ... file line.
**  Format assumed: 1:PRO_2:HA
**
**  After reading the file is **not** flushed to the next line.
**  Returns status flag.
*/

static Status ReadAtomData_biosym( FILE *fp, char *atnam, char *resnam, 
                            int *atnum, int *resnum, const int size )
{
    int   c;
    int   count;
    
    while ( ( c = getc( fp ) ) > 0 )
        if ( c == ':' )                 /* skip to residue name */
            break;
    count = size;
    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( c != '_' )
        {
            if ( !(--count) )
            {
                puts( "Exceeding allocated size of residue name string" );
                return( read_err );
            }
            *resnam++ = c;              /* store residue name */
        }
        else
            break;
    }
    if ( fscanf( fp, "%i", resnum ) < 1 )
    {
        puts( "Error in reading residue number from Biosym file" );
        return( read_err );
    }                                   /* read residue number */
    while ( ( c = getc( fp ) ) > 0 )
        if ( c == ':' )                 /* skip to atom name */
            break;
    count = size;
    while ( ( c = getc( fp ) ) > 0 )
    {
        if ( !IsEmptyChar( c ) )
        {
            if ( !(--count) )
            {
                puts( "Exceeding allocated size of atom name string" );
                return( read_err );
            }
            *atnam++ = c;               /* store atom name */
        }
        else
            break;
    }
    *atnum = 0;

    return( status_ok );
}

/****************  ReadCoordData_biosym  ****************/
/*
**  Read coordinate information from one Biosym car file line.
**  Format assumed: free format floats.
**
**  After reading the file is **not** flushed to the next line.
**  Returns status flag.
*/

static Status ReadCoordData_biosym( FILE *fp, float *x, float *y, float *z )
{
    if ( fscanf( fp, " %f %f %f", x, y, z ) < 3 )
    {
        puts( " * Error in reading coordinate data from Biosym file" );
        return( read_err );
    }
    return( status_ok ); 
}

/****************  MDFStorePseudo  ****************/
/*
**  Store a pseudo atom description from a Biosym mdf file.
**  Format assumed:
**  @list pseudoatom  PCOL:PRO_4:HDX A
**
**  alternative code: scanf("[^:]",...)
*/

static void MDFStorePseudo( char *line, Psatom *psatom )
{
    char *c;
    
    while ( *line++ != ':' )
        ;                               /* skip to residue name */
    c = psatom->res_nam;
    while ( *line != '_' )              /* store residue name */
        *c++ = *line++;
    if ( sscanf( ++line, "%i", &(psatom->res_num_orig) ) < 1 )
    {
        printf( "\n * Error in reading mdf file\n" );
        exit( 1 );
    }                                   /* read residue number */
    while ( *line++ != ':' )
        ;                               /* skip to pseudo name */
    c = psatom->pseudo_nam;
    while ( !IsEmptyChar( *line ) )     /* store pseudo name */
        *c++ = *line++;
    
    PRINTID( psatom->res_nam );
    printf( " %i ", psatom->res_num_orig );
    PRINTID( psatom->pseudo_nam );
    putchar( '\n' );

}

/****************  MDFStoreAtom  ****************/
/*
**  Store one (or more: recursive call) constituent atoms from the line(s).
**  following a pseudo atom description in a Biosym mdf file.
**  Format assumed:
**  PRO_4:HD1 PRO_4:HD2 
*/

static void MDFStoreAtom( char *line, Psatom *psatom )
{
    int   i;
    char  *c;
    
    if ( ( i = (psatom->num_form)++ ) == max_PSFORM )
    {
        printf( "\n * Out of array bound, change max_PSFORM\n" );
        exit( 1 );
    }
    while ( *line++ != ':' )
        ;                               /* skip to atom name */
    c = psatom->form_nam[i];
    while ( !IsEmptyChar( *line ) )     /* store atom name */
        *c++ = *line++;

    printf( " #%i:", i );
    PRINTID( psatom->form_nam[i] );

    if ( !IsEmptyStr( line ) )
        MDFStoreAtom( line, psatom );
}
